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ABSTRACT 

We consider a model of galaxy clusters in which the hot gas is in hydrostatic 
equilibrium and maintains energy balance between radiative cooling and heating 
by thermal conduction. We analyze the thermal stability of the gas using a 
Lagrangian perturbation analysis. For thermal conductivity at the level of ~ 
20 — 40% of Spitzer conductivity, consistent with previous estimates for cluster 
gas, we find that the growth rate of the most unstable global radial mode is ~ 6—9 
times lower than the growth rate of local isobaric modes at the cluster center in 
the absence of conduction. The growth time in typical clusters is ~ 2 — 5 Gyr, 
which is comparable to the time since the last major merger episode, when the 
gas was presumably well mixed. Thus, we suggest that thermal instability is not 
dynamically significant in clusters, provided there is an adequate level of thermal 
conduction. On the other hand, if the heating of the gas is not the result of 
thermal conduction or any other diffusive process such as turbulent mixing, then 
the thermal instability has a growth time under a Gyr in the central regions of the 
cluster and is a serious threat to equilibrium. We also analyze local nonradial 
modes and show that the Lagrangian technique leads to the same dispersion 
relation as the Eulerian approach, provided that clusters are initially in strict 
thermal equilibrium. Because cluster gas is convectively stable, nonradial modes 
always have a smaller growth rate than equivalent radial modes. 



Subject headings: galaxies: clusters — cooling flows — X-rays: galaxies — con- 
duction — hydrodynamics — instability 
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1. Introduction 

Galaxy clusters contain a large amount of hot diffuse gas radiating prolifically in ther- 
mal X-rays. In the absence of heat sources, the radiative cooling due to this emission should 
induce a subsonic inflow of gas in the central regions, leading to substantial condensation 
of cold gas in rich clusters (e.g., Fabian 1994 and references therein). However, recent high- 
resolution X-ray data from XMM-Newton and Chandra reveal no evidence for such cooling 
flows, nor any significant evidence for mass dropout. In particular, there is no evidence for 
gas at temperatures below about a third of the average temperature (Peterson et al. 2001, 
2003; Tamura et al. 2001; Bohringer et al. 2001; Fabian et al. 2001a; Molendi & Pizzolato 
2001; Matsushita et al. 2002; Johnstone et al. 2002), suggesting that there must be some 
heat source (or sources) that prevents the gas from cooling below this limit. Candidate heat- 
ing mechanisms include radiative and mechanical power from active galactic nuclei (AGN) 
(Ciotti & Ostriker 2001; Churazov et al. 2002; Briiggen & Kaiser 2002; Reynolds, Heinz, & 
Begelman 2002; Kaiser & Binney 2003), thermal conduction from the hotter outer regions of 
the cluster to the center (Narayan & Medvedev 2001; Fabian, Voigt, & Morris 2002; Gruzinov 
2002; Zakamska Sz Narayan 2003, hereafter ZN03), or perhaps both (e.g., ZN03, Ruszkowski 
& Begelman 2002; Brighenti & Mathews 2002, 2003). 

The effects of conduction in clusters have been widely discussed by many authors (e.g., 
Binney & Cowie 1981; Tucker & Rosner 1983; Bertschinger & Meiksin 1986; Bregman & 
David 1988; Gaetz 1989; Rosner & Tucker 1989; David et al. 1992; Pistinner & Shaviv 
1996; see Fabian 1994 and references therein). It was shown that thermal conduction, when 
unimpeded, can significantly reduce inferred mass deposition rates in cooling fiow clusters 
(Bertschinger & Meiksin 1986; Rosner & Tucker 1989), although such systems may end up 
being almost isothermal (Bregman & David 1988). Nevertheless, thermal conduction was 
considered to be unimportant in galaxy clusters as it was beheved that, in the presence of 
magnetic fields, the cross-field diffusion coefficient would be neghgibly small. Since magnetic 
fields are ubiquitous in clusters (e.g., CariUi & Taylor 2002), it appeared that the effective 
isotropic conduction coefficient k in the presence of magnetic fields would be smaller by 
orders of magnitude than the classical Spitzer (1962) value, Ksp, of an unmagnetized plasma. 

Narayan & Medvedev (2001) recently showed that cross-field diffusion of electrons is 
quite efficient if magnetic fields are fully turbulent and have a wide range of coherence length 
scales. Extending earher work by Rechester & Rosenbluth (1978), Chandran &: Cowley 

(1998), Chandran et al. (1999), and Malyshkin & Kulsrud (2001), and adopting the model 
of Goldreich & Sridhar (1995, 1997) for MHD turbulence, Narayan & Medvedev (2001) 
estimated that / = n/usp ~ 0.2 in a turbulent magnetized plasma. Very recently, Cho et 
al. (2003) have used direct numerical simulations to show that the turbulent diffusion of a 
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scalar field in MHD turbulence is as efficient as, or perhaps even more efficient than, the 
prediction of Narayan & Medvedev (2001). These studies have led to the revival of the idea 
that conduction might be an important heating source in clusters. 

The viability of thermal conduction as a heating mechanism has found support in the 
work of Voigt et al. (2002), Fabian, Voigt, & Morris (2002), and ZN03, who showed that the 
level of condiictivity needed to fit the observed temperature distributions of X-ray clusters 
is consistent with the theoretical estimate of Narayan & Medvedev (2001). In particular, 
ZN03 explicitly solved the equations for hydrostatic equilibrium and energy balance between 
radiative cooling and conductive heating and showed that conduction with / ~ 0.2 — 0.4 
provides reasonable profiles of density and temperature for several clusters. A similar result 
was also obtained by Brighenti & Mathews (2003) using numerical simulations. Although 
other heating sources, including AGN jets, cannot be ruled out, these studies suggest that 
an energy balance based purely on conduction is not unreasonable. 

While the above studies show that it is possible to construct equilibrium cluster models 
with conductive and/or AGN heating, there is no guarantee that the equilibrium will be 
stable. This is because hot, optically-thin. X-ray-emitting gas is well-known to be thermally 
unstable (Field 1965). Since the growth time of the thermal instabihty is comparable to the 
cooling time (e.g., Fabian 1994), one might expect rapid mass dropout as a result of the 
instability, even when there is a source of heat to eliminate the classic cooling flow. On the 
other hand, the absence of any evidence for gas below a few kcV in clusters implies that the 
instability is either absent or very slow. A likely reason for the lack of instability is thermal 
conduction, which is known to suppress thermal instability on small scales (see, e.g.. Field 
1965; Defouw 1970; Malagoh et al. 1987; McKee & Begelman 1990; ZN03). However, all 
previous analyses of this process have been limited to local WKB-type perturbations, where 
the wavelength of the perturbations is much smaller than the local radius of the system. ZN03 
found from their WKB analysis the intriguing result that perturbations with wavelengths 
up to almost the radius are stable, but that longer wavelengths are probably unstable. Such 
large-scale variations can be analyzed only through a full global mode analysis, which has 
not been done so far. 

In this paper wc present a detailed and formal analysis of local and global modes of 
thermal instability in equilibrium galaxy clusters with conduction. In §2, we describe the 
basic equations we solve and present equilibrium solutions for our fiducial model cluster. 
In §3, we carry out a linear stability analysis of the assumed equilibrium and calculate the 
global radial modes of the system. We show that there is a single unstable mode, whose 
growth rate is much lower than the usual thermal instability growth rate in the absence of 
conduction. In §4, we confirm the main results by means of numerical simulations which show 
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the development of the thermal instability under various assumptions. In §5, we consider 
the stability of nonradial modes and clarify a few points on which there has been confusion 
in the literature. We conclude in §6 with a brief summary of the results and a discussion of 
the implications. 



2. Perturbed Equations 

2.1. Basic Equations 

We consider the thermodynamic evolution of a hot plasma subject to radiative cooling 
and thermal conduction. We do not consider any dynamical effects of magnetic fields. We 
also neglect the self-gravity of the gas which is generally weaker than the gravity of the dark 
matter.^ The governing hydrodynamic equations are 

| + pV-v = 0, (1) 

^ + -VP + V$ = 0, (2) 
at p 

1 ''f T ^*+,£ + v.F = 0. (3) 



— 1 dt ^ — 1 p dt 

Here, d/dt = d/dt + v • V is the Lagrangian time derivative, p is the mass density, v is 
the velocity, T is the temperature, $ is the gravitational potential, 7 = 5/3 is the adiabatic 
index of the gas, pC = 2.1 x 10~^'^ngT^/^ erg cm~'^ s~^ is the energy loss rate per unit volume 
due to thermal bremsstrahlung (Rybicki & Lightman 1979; ZN03), and P is the thermal 
pressure for which we adopt the equation of state of an ideal gas, 

P = = —UekBT, (4) 

pruu p 

where m„ is the atomic mass unit, rie is the electron number density, and p and pe denote the 
mean molecular weight per hydrogen atom and per electron, respectively. We use p — 0.62 
and Pe = 1.18, corresponding to a fully ionized gas with hydrogen fraction X — 0.7 and 
helium fraction Y = 0.28 (ZN03). 

In equation (3), the conductive heat flux F is determined by 

F = -kVT, (5) 



^We found that the exphcit inclusion of self-gravity increases the growth rates of global thermal instability 
by only ~ 3 % from the non-self-gravitating values. 
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where the conductivity k is a fraction / of the classical Spitzer (1962) conductivity Ksp, 

« = /«Sp = / ^— r ergs K cm , (6) 

in7\.c 

with the Coulomb logarithm In Ac ~ 37. We do not consider nonlocal heat transport as in 
the model of Chun & Rosner (1993). In an unmagnetized plasma, conduction is isotropic 
and / is unity. The presence of a magnetic field generally reduces / by resisting motions of 
thermal elections across the field lines. Although the suppression of k would be very high 
if magnetic fields are uniform or only moderately tangled (e.g., Chandran & Cowley 1998), 
Narayan & Medvedev (2001) showed that / ~ 0.2 in a fully turbulent plasma medium in 
which magnetic fields are chaotic over a wide range of length scales. In this paper, we assume 
that / is constant in both space and time. 

We use the NFW form of the dark matter distribution p^M to determine the gravitational 
potential $ through 

= AnGp^r) = (7) 

(r + rc)(r + r^)^ 

with a softened core radius Tc, a scale radius r^, and a characteristic mass Afo- Utilizing the 
mass-temperature relation of Afshordi & Cen (2002) and the mass-scale relation of Maoz 
et al. (1997), we determine Mq and for a given cluster from the observed temperature 
in the outer regions of the cluster (see ZN03). For Tc, which determines the shape of the 
potential in the very inner parts, we adopt the best-fit values (either Tc — Q or Tc — rs/20) 
recommended by ZN03. 

Without any heating source to compensate for X-ray cooling, a cluster would lose its 
thermal energy at a rate pC Assuming that cooling occurs at fixed pressure, we may define 
from equation (3) the isobaric cooling time^ 

Similarly, we define the conduction time as 



^Note that the cooling time as defined here is yUe/A* times longer than the conventional definition, 
5fcBT/2neA, with A = pC/n^ (of. Sarazin 1988; David et al. 2001). At constant pressure, equation (3) 
gives T{t)/T{0) = ne(0)/ne(t) = (1 - 3t/2Uoi)'/'. 
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where A is the length scale over which the temperature changes appreciably. By comparing 
equations (8) and (9), we see that conductive heating becomes comparable to radiative 
cooling at A ~ Ai;' where 

is the Field length (Field 1965; McKce & Bcgclman 1990). This corresponds to the length 
scale below which thermal conduction erases temperature perturbations completely (Field 
1965). Since the length scale for temperature variations is typically ~ 10 — 100 kpc, we con- 
clude that thermal conduction with / ^ 0.1 can potentially support galaxy clusters against 
radiative coohng (see Narayan & Medvedev 2001; Fabian, Voigt, & Morris 2002; Gruzinov 
2002; ZN03). 



2.2. Initial Equilibrium 

It is straightforward to find initial equilibrium solutions of equations (l)-(7). We assume 
that the equilibrium is spherically symmetric, static, and time-independent. Working in 
spherical coordinates (r, 6, (f)), we simplify equations (2), (3), and (5) to 

— = -p— (11) 

dr dr ' 

~ ir'P^) - -PC (12) 



where Fr denotes the radial heat flux. 



Full numerical solutions of the differential equations (7) and (11)-(13) were presented by 
ZN03. They varied / to find the solutions that give the best fits to the observed density and 
temperature profiles of ten clusters. They found that thermal conduction with / ~ 0.2 — 0.4 
explains the density and temperature distributions of five clusters fairly well (A1795, A1835, 
A2199, A2390, RXJ1347.5-1145). However, five other clusters require unphysically large 
values of /, indicating that those clusters are incompatible with a pure conduction model 
and require other heat sources. In this paper we take the best-fit parameters obtained by 
ZN03 for the five clusters that are consistent with conduction, and we analyze their stability 
to global radial modes. Table 1 lists the model parameters and various time scales. We 
adopt the cluster A1795 as our fiducial model. 



- 7- 



Figure 1 illustrates an equilibrium solution for A1795. The temperature is minimum at 

r = and begins to rise sharply at r ~ 5 kpc. Expansion of the variables near the center 
gives T(r) ~ Tq + {pqCq/Quqit'^ , where the subscript "0" indicates the values at the center. 
Therefore, stronger cooling or smaller conductivity would cause the temperature to increase 
faster. Notice that the temperature is always an increasing function of the radius. Radial 
heat influx is largest at r ~ 50 kpc beyond which the low density makes cooling as well as 
conductive heating unimportant. Figure Id plots the local isobaric cooling time ^cooi defined 
in equation (8) (solid fine) and the growth time, (see eq. [27]), of local isobaric thermal 
perturbations without conduction (dashed line), which are compared with the growth time 
of the global radial mode in the presence of conduction (igrow; dotted line): we discuss these 
instability time scales in the next section. 

The radial equilibrium profiles of A1795 presented in Figure 1 are close to but not exactly 
the same as those derived by Ettori et al. (2002). Since A1795 is known to contain a cool fila- 
ment near the center (Fabian ct al. 2001b), the real density and temperature distributions are 
somewhat non-axisymmctric. Our axisymmetric equilibrium model of A1795 should, there- 
fore, be regarded as an idealized version of what is a much more complex distribution of the 
intracluster medium in A1795. As we will show below, however, the lowest-order radial mode 
of thermal instability is the fastest growing mode, so that neglecting the non-axisymmetric 
parts in the background profiles probably does not affect the results significantly. 



3. Linear Analyses of Radial Modes 

3.1. Lagrangian Perturbations 

We linearize equations (l)-(5) in the Lagrangian framework. A Lagrangian perturbation, 
represented by an operator A, is related to an Eulerian perturbation 5 in the usual way, 

A = (5 + e-V, (14) 

— * 

where the vector ^ measures the Lagrangian displacement of a fiuid element from its unper- 
turbed location (e.g., Shapiro & Teukolsky 1983). One of the advantages of adopting the 
Lagrangian approach is that it simplifies the perturbed equations greatly. For instance, the 
perturbed continuity equation (1) becomes 

Ap=-pV-e (15) 

Nevertheless, both Lagrangian and Eulerian descriptions should give the same results, espe- 
cially if the initial state is static and in complete equilibrium. Various properties of A and 
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the commutation relations associated with A and 5 can be found in Shapiro & Teukolsky 
(1983). 

We assume that the perturbations are all radial. Applying A to equations (2)- (5) and 
using equation (15), we obtain 



d% P d f 1 dr^i. 
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(17) 
(18) 



where Ct = dC/dT\p, Cp = dC/dp\T-, is the radial component of ^, and Lr = —iirr'^Fr is 
the radial heat luminosity. We use and AL,. as independent variables. 

We seek hnear eigenmodes that behave as ~ e"'* with time. We then rewrite equations 
(16)-(18) as 
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(19) 

(20) 

(21) 



which are our desired perturbation equations. This set of ordinary differential equations 
forms an eigenvalue problem for global modes, with a as the eigenvalue, which can be solved 
numerically subject to appropriate boundary conditions. Before considering the full problem, 
we study local modes that do not depend on boundaries to gain some physical insight. 



3.2. Local Radial Modes 

Let us consider local WKB perturbations of the form ~ g^fcrr+crt ^.^^ jg^. assume that 
the radial wavenumber kj. satisfies hf-r ^ 1, kr{d\nP/dr) ^ 1. We may then neglect the 
spherical geometry and ignore local gradients of the background quantities relative to spatial 
gradients of the perturbed variables. We also assume that a^ < Pkj/p, corresponding to 
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slowly evolving perturbations; this eliminates sound waves from consideration. In this local 
approximation, equation (19) simplifies to 

AT 

while equations (20) and (21) may be combined to give 

-{Pa - p'Lp - pC)iKir = + P^'^-r + 

Eliminating and IST /T from equations (22) and (23), we find 



AT 



(23) 



<J^<Joo-^^kl (24) 



where (Too, defined by 



-l-^(2r^r Tr \ l-lpT(dC/T\ 

is the growth rate of isobaric thermal perturbations without conduction (Field 1965). 

Equation (24) is a local dispersion relation for thermal fluctuations.'^ In the absence 
of conduction, a = a^o, so that instability develops if the cooling function satisfies the 
generalized Field criterion 

[dC/T\ 



\ dT ) 



< 0, (26) 



p 



for isobaric thermal instability (Balbus 1986). For X-ray emitting clusters of galaxies with 
jC oc pT^f"^, the condition (26) is easily met. The corresponding growth time amounts to 

-1 2 

^OO = = g^COob (27) 

suggesting that local radial disturbances will grow slightly faster than the isobaric cooling 
time of the system. Figure Id plots too for A1795 as a dashed fine. 

It is well known that thermal conduction stabilizes short-wavelength perturbations 
against thermal instability (e.g.. Field 1965; Malagoli et al. 1987; McKee & Begelman 
1990). It not only reduces growth rates but also suppresses thermal instability completely 



^Equation (24) can be obtained by taking the limit of kr » po^ jP in the dispersion relation (15) of Field 
(1965) for a uniform medium. Note, however, that the condition £ = in his initial equilibrium state gives 
(Too = —[(7 — l)p/7-P](5£/9T)p instead of equation (25). 
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if kr > (3/2)^/^/cir, where kp — 27:/ Xp- The existence of unstable modes, therefore, requires 
large wavelength perturbations. ZN03 showed that the clusters they studied are thermally 
stable for all wavelengths up to approximately the radius. The fate of modes with longer 
wavelengths was unclear in their study since a local analysis is no longer valid. A global 
analysis is required to study such large scale perturbations. This is the topic of the next 
subsection. 



3.3. Global Solutions 

To analyze the global stability problem, we solve equations (19)-(21) numerically. Since 
these equations are equivalent to four first-order differential equations, we need to specify 
four boundary conditions. The two inner boundary conditions are 

^ = and ^ = 0, atr = 0. (28) 
dr 4:TTr^ 

The first condition guarantees that the solutions are regular, while the second condition 
corresponds to a zero gradient in the perturbed temperature at the center. Since = at 
r = for radial modes, we fix ^r/i^ = 1 as a normalization condition. 

The two other boundary conditions come from the outer boundary which we arbitrarily 
locate at Tf, = 1 Mpc. Since all cigcnfunctions that we have obtained are found to decay 
rapidly with increasing radius, the solutions are quite insensitive to the particular choice of 
the outer boundary conditions as well as the location of the outer boundary. This is because 
the cooling time well exceeds the Hubble time in the outer parts of clusters (Fig. Id), and so 
the thermal instability develops very slowly near the outer boundary. The results presented 
in this paper correspond to the specific boundary conditions 

= and AT = 0, at r = n. (29) 

We have tried various other outer boundary conditions such as a fixed pressure, a fixed 
entropy, etc. The results are unchanged. 

The numerical calculations proceed as follows. We first fix a and integrate equations 
(19)-(21) from r = to r = r;,, setting a = AT/T\r=o to an arbitrary value. At the outer 
boundary, we check the first condition in equation (29), update a using the Newton- Raphson 
method, and continue iterating until convergence is attained. Then, we scan a in the range 
10~^ < \a\rb/co < 10^, where cq is the adiabatic sound speed at the cluster center, and use 
the second condition in equation (29) as a discriminant for solutions. By this procedure, we 
do not miss any eigensolution. We allow for both real and complex a. 
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When a is real and positive, we find that a cluster with / ^ 0.2 has only one unstable 
global mode. In Figure 2, we plot the eigenfunctions of the unstable mode for all five clusters 
listed in Table 1. The perturbations have largest amplitude near the center (r < 1 kpc) and 
decay rapidly with increasing radius. This is consistent with the notion that the central 
regions have a much shorter cooling time and are thus more prone to thermal instability. 
Note that ^^/^ has the same sign as AT/T, indicating that the thermal instability leads to 
a mass inflow. We will show in §4 that the mass inflow rate driven by the instability is, 
however, much lower than that in conventional cooling flows. As the last two columns of 
Table 1 shows, the growth time of the most unstable global radial modes is typically ~ 2 — 5 
Gyr, which is about 6 to 9 times longer than that of local isobaric modes at the cluster center 
in the absence of conduction. This time scale is not much shorter than the age of the system 
since the last major merger. Therefore, thermal instability in the presence of conduction is 
unlikely to be a serious threat to clusters. 

Why is the growth time tgrow of the mode so much longer than the growth time too,o of 
the local isobaric mode at the cluster center? The reason is that the eigenfunction extends 
over a considerable range of r, of order tens of kpc (Fig. 2). Over this range, too increases 
from too,o the center to much larger values (Fig. Irf). Since the growth time of the mode 
corresponds to an appropriate eigenfunction- weighted volume-average of too, and since the 
averaging is dominated by r ~ tens of kpc, this causes tgrow to be significantly larger than 

^00,0- 

How cfTcctivcly docs conduction stabilize thermal instability? To address this question, 
we have solved equations (19)-(21) using different values of / for the perturbations. We adopt 
the same density and temperature distributions as shown in Figure 1, which means that we 
take the same / = 0.2 for the initial unperturbed state. However, we allow / to vary in the 
perturbations. In Figure 3, we show the resulting mode growth rates a relative to (Joo^o, the 
central value of the growth rate of local thermal instability without conduction (see eq. [25]). 
When / > 5 (although unphysical), we find that conduction quenches thermal instability 
completely; that is, there is no unstable mode. For 0.14 </< 5, there exists a single 
unstable mode that has no node in the eigenfunction for t^rl'i'- The case shown in Figure 2, 
/ = 0.2, corresponds to this range. When / is below 0.14, a new mode appears and now the 
system has two unstable modes. Because the new mode has a node in ^^jr somewhere inside 
the cluster, its growth rate is smaller than the zero-node mode, as equation (24) suggests. 
As / keeps decreasing, new modes having lower growth rates successively emerge. In the 
limit of zero conduction, there is an infinite number of unstable modes, with the maximum 
growth rate equal to a"oo,o- 

Figure 3 shows that over much of the range of /, the growth rate increases very slowly as 
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/ decreases. As before, the reason may be traced to the fact that the eigenfunction extends 
over a range of r and a represents a suitable average of a^o over this volume. Thus, even for 
/ = 0.01, the fundamental n = mode has a growth rate about 4 times smaller than aoo,o- 
Eigenf unctions become more centrally concentrated as / decreases, and thus a mode with 
smaller / generally has a larger growth rate, though the increase of a with decreasing / is 
rather slow (Fig. 3). 

In addition to unstable modes for which thermal conduction plays a stabilizing role, 
clusters also contain decaying modes with real and negative a. Figure 4 displays examples 
of cigcnfunctions for a few selected decaying modes and shows the dependence of a on the 
reduction factor /. In contrast to unstable modes, which have smaller frequencies as the 
number of nodes increases, higher-order decaying modes have larger (negative) frequencies. 
This is consistent with equation (24)^, which implies that the presence of decaying modes is 
a simple manifestation of diffusion ironing out temperature perturbations. 

We have searched for global solutions with complex cr, i.e., overstable modes, but found 
no overstable radial mode with Re((T) > under the imposed boundary conditions. Phys- 
ically, thermal overstability occurs when sound waves are amplified by absorbing (losing) 
heat during the compression (rarefaction) phase of oscillations (Field 1965). Since the cool- 
ing function for X-ray emitting clusters does not satisfy the isentropic instability criterion 
of Field (1965), the absence of radial overstability is in fact guaranteed. This leaves pure 
unstable modes as the sole growing modes. Although some undcrdamping modes (complex 
a with negative real parts) were found, they are of no interest and we do not discuss them 
further. 



4. Nonlinear Evolution of Radial Modes 

As a check of the global linear stability analyses described in the previous section, we 
have solved the full dynamical equations (l)-(7) using a time-dependent approach. Starting 
with the equilibrium density and temperature profiles of A1795 calculated in §2.2, we have 
run three models: (model A) a cluster without any heating; (model B) a cluster with no 
conduction, but constantly heated (by some fictitious agency) at a rate such as to maintain 
thermal balance in the initial equilibrium; (model C) a cluster with self-consistent conductive 
heating with / = 0.2. Model A simulates the standard cluster coohng fiow problem, whereas 
models B and C are intended to explore the growth of thermal instability from an initial 



"^Since a cx —nTk^ for decaying modes, equation (24), which assumes <C Pk^/p, does not describe 
higher-order modes very accurately. 
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equilibrium state, with and without thermal conduction. Note that our adopted coohng 
function for X-ray free-free emission becomes invahd when the gas in models A-C cools 
below ~ 1 keV. Nevertheless, these simulations allow us to confirm the growth rates of 
thermal instability that we computed in §3 and to estimate mass inflow rates resulting from 
either radiative cooling or thermal instability. 

We follow the nonlinear evolution of om model clusters using the ZEUS hydrodynamic 
code (Stone & Norman 1992). We construct a logarithmically-spaced radial grid, with 500 
zones, from r = 1 kpc to 1 Mpc, and carry out a one-dimensional simulation. We implement 
a fully implicit algorithm for thermal conduction in the energy equation (Press ct al. 1992). 
All the variables at both the inner and outer boundaries are assigned to have a vanishing 
gradient across the boundaries, except for the radial velocity. We fix the radial velocity to 
be zero at the outer boundary, while allowing it to vary as a linear function of radius at 
the inner boundary. On the initial equilibrium, we add as a perturbation the most unstable 
global eigenfunction for density that we found in §3.3, with an amplitude equal to 0.1% of 
the background density at the cluster center. 

Figure 5 shows the evolutionary histories of the maximum density and the mass inflow 
rate M = —Anr^pVr at r = 10 kpc for models A-C. The two dotted hues in Figure 5a 
represent the analytic estimates for the growth rates of thermal instability with and without 
conduction; the predictions arc in excellent agreement with the results of the numerical 
simulations. We have also compared the perturbed density and velocity profiles in the linear 
regime of model C with the corresponding eigenfunctions obtained from the linear theory. 
Even though the radial grid does not cover the inner 1 kpc, they agree to within ~ 3%. Radial 
profiles of electron number density, radial velocity, and temperature for a few selected epochs 
are plotted in Figure 6. 

Model A, which is not in equilibrium and is subject to strong radiative cooling, immedi- 
ately experiences radial infiow of mass everywhere. Although the temperature of the cluster 
decreases steadily as the gas cools, the thermal pressure does not drop owing to adiabatic 
compression, which in turn maintains the radial velocity at the level of about ~ 1 — 2% 
of the central sound speed until t ~ 0.6 Gyr. During this period, we find Ue oc T~^-°^ so 
that the cooling flow is nearly isobaric. As the insert in Figure 5a shows, the maximum 
density follows the prediction of isobaric cooling fairly well. When the cooling modifies the 
density and temperature distributions significantly {t > 0.6 Gyr), the increased cooling rate 
in the central parts induces larger mass infiows, causing the central density to increase in a 
runaway fashion. The change of density and temperature is so rapid that there is no chance 
for thermal instabihty to grow. We stop the simulation of model A ai t — 0.85 Gyr when 
the central temperature becomes vanishingly small. The mass infiow rate depends rather 
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sensitively on time and is a linear function of radius for r < 250 kpc, with a slope of ~ 30 Mq 
yr"^ kpc~^, at the end of the run (see Figure 56). 

Models B and C are initially in equilibrium with heating balancing cooling, but they 
both experience thermal instability. Without conduction, perturbations in model B grow 
at a rate of (7oo,o = (0.64 Gyr)~^ and become fully unstable within ~ 4 Gyr. Although the 
growth time of thermal instability in model B is slighter shorter than the isobaric cooling 
time in model A (see eq. [27]), model B takes a longer time for eventual runaway than model 
A. This is simply because the former requires the growth of perturbations from a small 
amplitude, while the latter changes its background state immediately. 

Model C clearly illustrates the stabilizing effect of thermal conduction. With / = 0.2, 
conduction causes the maximum growth rate of thermal instability to be about six times 
smaller than in the nonconducting model B, thereby delaying runaway growth until t ~ 
24 Gyr (for the particular amplitude of initial perturbations). The mass inflow rate is 
correspondingly very small throughout its entire evolution. Had we started the evolution 
with a higher initial perturbation amplitude, it would of course take less time to reach a fully 
nonlinear state. Small-scale perturbations in real clusters may well have high amplitudes, 
but they are readily erased by conduction. Since it seems unlikely that the amplitudes 
of perturbations that arc spherically symmetric and coherent over 100 kpc are more than 
~ 10% of background quantities, we believe that galaxy clusters with / ^ 0.2 are not likely 
to exhibit a strong thermal instability for many Gyr. 

5. Local Nonradial Modes 

We now carry out a local analysis of nonradial Lagrangian perturbations in a static, 
stratified medium. Wc do not attempt to solve the full eigenvalue problem for nonradial 
modes, but instead focus on the local behavior of nonradial modes in the presence of den- 
sity and/or temperature stratification. Similar work has been reported by White & Sarazin 
(1987) and Malagoli et al. (1987) who used Eulerian perturbations and showed that a restor- 
ing buoyancy force from a stable entropy gradient can change thermal instability into an 
overstability, confirming the previous results of Defouw (1970). Using Lagrangian perturba- 
tions, on the other hand, Balbus (1988) showed that the negative radial gradient of a net 
cooling function leads to convectively unstable flows in a non-equilibrium system. He further 
claimed that the Eulerian and Lagrangian approaches give different results, although in a 
subsequent paper Balbus & Soker (1989) showed that the Lagrangian approach of thermal 
instability in dynamical media produces results that are apparently consistent with those 
from the Eulerian analysis (Balbus 2003, private communication). We show below that the 
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Lagrangian approach in fact gives exactly the same dispersion relation for nonradial thermal 
instability as the Eulerian analyses referred to above, provided the background is in initial 
thermal equilibrium. 

We use the Lagrangian technique and linearize equations (l)-(5) assuming that pertur- 
bations have small amplitudes of the form ~ Yim(^, 0)e''^''''+'^*, where Yijn{0, 0) is the usual 
spherical harmonic. The detailed procedure to derive a local dispersion relation is given in 
the Appendix: we write only the final result here, 

' - - V^* ) " + s^"-" ^ ^rm^"" + ^ ■ = 0, (30) 

where kt = [/(/ + l)]^/^/r is the tangential wavenumber, k'^ = k'^ + and 

2 _ d$ fldlnP d\np\ 

= Iv7^^^ d^J ^ ^ 

is the Brunt- Vasaila frequency for convective motions in a radially stratified atmosphere. 
When the conductivity k = 0, equation (30) reduces to equation (3.17) of Balbus (1988)^. 
If d{pC + 'V • F)/dr = in the last term, equation (30) becomes exactly the same as the 
dispersion relation of the Eulerian approach (Malagoli et al. 1987). 

It is the last term in equation (30) that has led to the discrepancy between the Eulerian 
and Lagrangian perturbation analyses. This term vanishes if the background is in thermal 
balance, either in the presence of thermal conduction in which case pC + V ■ F = 0, or with 
some other sources of heating in which case C = everywhere. If one perturbs equation 
(3) using the Eulerian operator 6 and Lagrangian operator A, respectively, and takes the 
difference of the resulting equations, the residual term is ^ ■ V{pC + V ■ F). This explains 
the origin of the last term in equation (30) and the apparent discrepancy of the Eulerian 
and Lagrangian approaches. The two approaches are consistent with each other only if the 
background state is in strict thermal equilibrium. 

One may argue that the assumption of initial thermal equilibrium is very special and 
a local analysis can be carried out even in a non-equilibrium system. However, such an 
analysis would be misleading unless the change of the background state either is allowed 
for in doing the perturbations or occurs slowly compared to the growth of perturbations 



^Since Balbus (1988) took pC^Q m an initial state, it is unclear how to compute an effective adiabatic 
index Fi in his notation. Applying A directly to his equation (2.1c), one obtains Fi = (5(t/2+t0,-.p)/(3(t/2+ 
Q + 6r,p) instead of his equation (2.8). Plugging this expression into his equation (3.17) would produce 
exactly the same local dispersion relation as our equation (30) when k = 0. Then, the condition dQ/dr < 
in his equation (3.28) is replaced by d{pC)/dr < 0. 



16 



of interest. As we showed in §3, the coohng time in galaxy clusters is comparable to the 
growth time of thermal instability. Moreover, in §4, wc performed a numerical simulation of 
a cluster that is initially out of thermal equilibrium (model A) and showed that the resulting 
cooling flow quickly dominates the cluster evolution, leaving no time for the development of 
a local thermal instability. 

Since we assume the system is in initial thermal balance, we rewrite equation (30) as 



which was already derived and discussed by Malagoli ct al. (1987). Compared to equation 
(24) for radial perturbations, equation (32) has an extra term that arises from the local 
buoyancy force. Since dT/dr > and dP/dr < 0, therefore u^y > 0, and galaxy clusters are 
stable to convective instability. The buoyancy term has then a stabilizing effect, reducing 
the growth rate of the thermal instability and even altering the character of the instability to 
an overstability when kt/k is high enough (Defouw 1970; Malagoli et al. 1987). If a medium 
is thermally overstable, a bubble that is displaced radially outward (inward) becomes more 
(less) dense compared to the adiabatic case, thus experiencing an increase of its oscillation 
amplitude with time. 

Prom equations (24) and (32), it is trivial to show that for a given wave- vector k, Re{a) 
of the local nonradial mode is always smaller than that of the local radial mode (e.g.. White & 
Sarazin 1987), implying that the latter is more unstable. This suggests that global nonradial 
modes, if they exist, would have smaller growth rates compared to global radial modes. 



Narayan & Medvedev (2001) and Cho et al. (2003) have shown that thermal conduction 
is quite efficient in a fully turbulent, magnetized plasma. The diffusion coefficient in such 
a medium is approximately a fraction / ~ a few tenths of the Spitzer coefficient (the full 
Spitzer coefficient applies to an unmagnetized medium) . Conduction at this level is sufficient 
to provide significant heat to the cooling gas in a galaxy cluster and to maintain energy 
balance in the cluster center (Narayan & Medvedev 2001; Voigt et al. 2002; Fabian, Voigt, & 
Morris 2002; ZN03). This is an attractive explanation for the suppression of mass dropout in 
cooling cores of clusters. However, the demonstration that one can build equilibrium models 
with heating equal to coohng is not sufficient, since the hot gas is hkely to be subject to a 
thermal instability (Field 1965). To be fully consistent with the new X-ray data on clusters, 
which show the absence of gas below a few keV, it is necessary to demonstrate that the 
thermal instability is either fully suppressed or is at least much weaker than expected. 




(32) 



6. Summary and Discussion 
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It is well-known that thermal conduction suppresses thermal instability in a hot gas for 
perturbations with short wavelengths. Indeed, using a local WKB-like analysis, ZN03 have 
shown that perturbations with wavelengths up to a fraction of the radius are thermally stable 
in clusters. However, this still leaves open the possibility that one or more long-wavelength 
global modes may be unstable. In this paper, we have analyzed in some detail the stability 
of such global modes. 

Following ZN03, wc set up the equilibrium density and temperature profiles of clusters, 
assuming strict hydrostatic and energy equilibrium (§2). By applying Lagrangian perturba- 
tions, wc derive a set of differential equations that describe the eigenvalue problem for global 
radial modes, with the mode growth rate a acting as the eigenvalue (eqs. [19]-[21], §3). In 
the absence of conduction (/ = 0), we find as expected that a cluster has an infinite number 
of unstable modes, all with rapid growth rates a ~ — (eqs. [24] and [27]), where t^o 
is the standard growth time scale associated with the thermal instability. For typical cooling 
flow clusters, the modes should grow in less than a Gyr, implying a serious threat to stability. 
However, when we set / = 0.2, roughly the value recommended by Narayan & Medvedev 
(2001) and ZN03, we find that all radial modes except one become stable. The one residual 
unstable mode is the fundamental nodeless (in ^r/''^) mode, i.e. the mode with the longest 
possible "wavelength" that can fit within the system. This mode is weakly unstable, with a 
growth time scale tgrow that is 6 to 9 times longer than too- Soker (2003) has recently posted 
a paper in which he describes a quasi-global stability analysis of clusters. However, both the 
initial equilibrium he assumes and the perturbations he considers are very approximate, and 
it is hard to compare his results with ours. 

For the clusters we have studied, tgrow for the lone unstable mode is ~ 2 — 5 Gyr (Table 
1). This is an interestingly long time scale, which is probably comparable to the elapsed 
time since the last major merger event in a hierarchical clustering scenario (typically ~ 7 
Gyr for massive clusters, e.g., Kitayama & Suto 1996). We imagine that major mergers (and 
perhaps to a lesser extent even minor mergers) leave the cluster gas in a highly mixed and 
turbulent state. If the merger drives the whole system well out of thermal equilibrium, mass 
dropout would occur very rapidly (< 1 Gyr) at the central parts (as in model A in §4), which 
may in turn allow the remaining gas to achieve a new thermal equilibrium that is similar to 
its present state. In a sense, the merger acts to reset initial conditions. After the merger, 
the unstable global radial mode would grow. However, if the growth time for the mode is 
comparable to the effective age of the cluster since the last merger, as appears to be the case 
for our models, then we do not expect the thermal instability to be a major problem. 

We have confirmed the results of §3 by means of numerical simulations of radial pertur- 
bations in a model cluster (§4, Fig. 5). Three models are discussed. Model A considers a 
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cluster without any heating to balance coohng. This model exhibits the classic coohng flow 
catastrophe, a violent and very rapid runaway at the cluster center. The runaway is so rapid 
that there is no possibility for the thermal instability to do anything. Model C considers 
a cluster in which cooling is exactly balanced by conductive heating in equilibrium. The 
numerical simulation shows that this model has a slow instability that grows at precisely the 
rate calculated in §3. For the particular initial conditions selected, the mode grows to the 
non-linear runaway stage only after 24 Gyr, which means that the model is for all practical 
purposes stable. Model B is an interesting in-between case in which heating and cooling 
are balanced, except that the heating is not from conduction but from some other local 
(unspecified) agency. This model shows the classic thermal instabihty. The perturbations 
grow quite rapidly, on the time scale too of the Field (1965) instability. In 4 Gyr, the model 
is completely destroyed by nonlinear runaway. 

We thus reach the following important conclusion from the simulations: apart from 
balancing heating and cooling, it is necessary also to make sure that the heating is of the 
right kind to control the thermal instability. Specifically, the heating must involve diffu- 
sive transport of energy so that short-wavelength perturbations are smoothed out and not 
allowed to grow, and even long-wavelength perturbations are partially stabilized. Thermal 
conduction is diffusive in nature and is quite effective in this respect, as we have shown in this 
paper. Turbulent mixing would be similarly effective. AGN heating via a jet might behave 
diffusively if the heat is transferred to the cluster gas via mechanical turbulence. Heating 
through mass infall in minor mergers might also be effective since in this case again the 
heat is hkely to spread through the gas via turbulence. However, radiative heating (Ciotti 
& Ostriker 2001) or cosmic ray heating (Loewenstein, Zweibel, & Begelman 1991) from an 
AGN are unlikely to control the thermal instability, since these heating mechanisms are not 
diffusive in nature. 

In this connection, we note that not all solutions to the cooling flow problem involve 
equilibrium models. Kaiser & Binney (2003) have described an interesting model involving 
AGN heating in which the gas goes through a limit cycle. For most of the time, the AGN 
at the center is in quiescence and the gas is not significantly heated. As the gas undergoes a 
cooling catastrophe, the gas density around the AGN increases, causing the AGN to switch 
on, to eject a powerful jet and to heat the cluster gas. The heated gas expands, the AGN 
switches off, and the gas starts the cycle again. The analysis we have carried out in this 
paper does not include such scenarios, since we have explicitly assumed that heating and 
coohng are balanced. 

In §5, we analyze the properties of nonradial perturbations in a cluster. This problem 
has been studied by Defouw (1970), White & Sarazin (1987), and Malagoh et al. (1987) 
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using an Eulerian approach, and by Balbus (1988) using a Lagrangian approach. We show- 
that the Lagrangian analysis leads to the same local dispersion relation (see eq. [32]) as the 
previous Eulerian analyses, provided that the background is in thermal equilibrium. For a 
cluster that is convectively stable (as all clusters are), u^y > (eq. [31]), and the buoyancy 
force associated with the entropy gradient plays a stabilizing role. As a result, if we consider 
radial perturbations with a given wave-vector that are thermally unstable, then non-radial 
perturbations with the same wave- vector always have a lower growth rate. This means that 
the most unstable mode is always a radial mode. Although this result is obtained via a local 
analysis, it is expected to be true for global modes as well. We have, therefore, not analyzed 
global non-radial modes. 

Is it an accident that the mode growth time that we have estimated is comparable to 
the age of the system as measured since the last major merger? We suggest that perhaps 
it is not. Imagine a cluster that is formed with sufficient gas such that several modes are 
initially thermally unstable with fairly short growth times. If the cluster were to start from 
a well-stirred initial state (say immediately after a merger), the various modes would grow 
and when the age of the system exceeds the growth time of any particular mode, that mode 
would go non-linear and cause a certain amount of mass to cool and drop out from the 
hot medium. The mass dropout will cause the amount of gas in the hot phase to decrease, 
and as a result, the other modes would become less unstable. With continued mass dropout, 
perhaps only one unstable mode would be left finally. This mode would remain in the system 
so long as its growth time is not much shorter than the current age. With increasing age, 
presumably this mode too will cause some mass dropout, but always in such a manner that 
the gas that is left will have a growth time comparable to the age. In this picture, mass 
dropout acts as a safety valve that enables the system to be always in a state of marginal 
equilibrium. If this scenario is valid, then it is of course natural that the clusters we observe 
have unstable modes with growth times of order several Gyr. 

In this paper, we have assumed the parameter / to be constant over the entire volume 
of a cluster. This may not be a good approximation in some clusters that exhibit sharp 
discontinuities in temperature and density (e.g., Markevitch ct al. 2000; Vikhlinin et al. 
2001a,b; Dupke & White 2003). These cold fronts, located at about r ~ 300 kpc, probably 
result from cluster mergers. Vikhlinin et al. (2001b) argued that strong magnetic fields 
parallel to cold fronts may be responsible for the low conductivity. If so, the existence of 
cold fronts is unhkely to affect the strong stabihzing role played by conduction for the bulk 
of the cluster gas. The cold fronts are presumably transient features, surviving for only a 
dynamical time (< 1 Gyr) before being disrupted by Kelvin- Helmholtz instability (Mazzotta 
et al. 2002) or by merger shocks (Nagai & Kravtsov 2003). 
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Prom spatial variations of temperature detected in the cluster A754, Markevitch et al. 
(2003) argued that, in addition to cold fronts, the bulk of the gas in this cluster may have 
conductivity much smaller than the Spitzer value. The critical question is how recently were 
the observed temperature inhomogeneities formed and how much longer will they last. If the 
cluster is undergoing a merger, for instance, the density and temperature fluctuations are 
all transient. Although conduction tries to smooth out local temperature fluctuations, these 
structures may be continuously created by subsequent mergers and/or the sloshing motion 
of the gas in the dark matter potential (Roettiger et al. 1998). Acoustic motions of the gas 
may also induce temporary fluctuations. 

Since a large fraction of galaxy clusters exhibit powerful extended radio sources at their 
centers, many authors have considered outbursts from the central AGN to be the source 
of heat to balance the radiative coohng in clusters. As noted above, this works best if the 
heating is done via some turbulent agency so that short wavelength thermal instability can 
be controlled. In addition, the AGN heating cannot be too strong since it would smooth out 
and erase the radial variations of mctallicity that have been observed in some clusters (e.g., 
Allen et al. 2001; Johnstone et al. 2002). 

A very real possibility is that both AGN heating and conductive heating work together 
in clusters, though their relative importance may vary from cluster to cluster (ZN03). It is 
even possible that AGN heating dominates in the inner regions, while conduction plays a 
more important role farther out. Recently, Ruszkowski & Begelman (2002) and Brighenti & 
Mathews (2003) demonstrated that simultaneous heating by AGN and thermal conduction 
produces quasi-static density and temperature profiles that are similar to those of observed 
clusters. In particular, Brighenti & Mathews (2003) showed that clusters with AGN heating 
alone exhibit either too high mass accretion rates or unrealistic temperature distributions, 
while conductive heating with / ~ 0.35 alone or together with AGN heating gives reasonable 
fits to observations. Although the effects of their computational methods (e.g., removal of 
cold gas out of the computational domain ) and prescriptions for the AGN heating are uncer- 
tain, their simulations suggest that the stabilization of the hot gas is most likely achieved by 
conduction. Therefore, even in clusters where AGN heating is dominant, the role of thermal 
conduction as a stabilizing agent should not be ignored (ZN03). 

We gratefully acknowledge helpful discussions with S. Balbus, L. David, G. Field, E. 
Ostriker, and N. Yoshida. We also thank an anonymous referee for useful comments. This 
work was supported in part by NASA grant NAG5-10780. 



A. Nonradial Modes 



In this appendix, we linearize equations (l)-(5) for nonradial modes of perturbations 
and derive a local dispersion relation for thermal instability in a stratified medium. As in §3 
for radial modes, our approach is Lagrangian and includes the effects of thermal conduction. 
Similar analyses have been performed by Malagoh et al. (1987) for Eulerian perturbations, 
and by Balbus (1988) for Lagrangian perturbations, but without conduction. The effects of 
a radial mass flow and magnetic flelds were studied by Balbus & Soker (1989) and Balbus 
(1991) using Lagrangian perturbations. We apply A to equations (2) and (3) and obtain 

= (V • + -WiPW • a - (^P^) - V(e • V$), (Al) 



J_|.AL, + «V? (at - - (V. . f.) V . F, (A2) 



where denotes the tangential component of ^ and V* is the tangential gradient 
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with 9 and (p denoting unit vectors in the polar and azimuthal directions, respectively. 
Equation (18) linking AT and AL^ is valid for both radial and nonradial perturbations. 

Following Balbus (1988), we assume that the perturbations are of the form (see also 
Cox 1980) 

e = Hr{r) + r6(r) VJF,™(^, 0)e'^*, (A4) 
where Yi^{9, (p) is the spherical harmonic, satisfying an identity 
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We assume the same angular and temporal dependences in the other perturbed variables 
AT and AL^; in what follows, we omit Yime"^* from all the perturbation variables. 



The tangential component of equation (Al) is integrated to give 
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Prom equations (A4) and (A6), we thus write 
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where the dimensionless parameter D, defined by 
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measures approximately the square of the ratio of the growth time to the sound crossing 
time across a tangential wavelength. 

Substituting equations (A4) and (A7) into the radial component of equation (Al), we 
obtain 
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Similarly, equation (A2) becomes 
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Note that when D = 1, corresponding to pure radial modes, equations (A9) and (AlO) are 
reduced to equations (19) and (20), respectively. Subject to proper boundary conditions, 
equations (A9), (AlO), and (21) may be integrated to yield solutions for global nonradial 
modes, but this is beyond the scope of the present paper. 

Let us define the local tangential wavcnumbcr kt = [/(/ + l)]^/^/r and the local radial 
wavenumber kr = dlnx/dr, where x refers to any perturbed variable. Let us also focus 
on local modes which vary rapidly in both radial and tangential directions and grow slowly 
compared to the sound crossing time across their wavelengths, i.e. krr ^ 1, kr{d\n P/dr) 3> 
1, and P^^ jP (or P> ^ 1)- Using equation (21) in the following alternative form. 



^ = |:(-Ar + F.e.)-e.(v-F), 



(All) 
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we simplify equations (A9) and (AlO) to 

dlnPdlnT 



k'-D- 



dr 



d 



dr 

dlnP 



dlnP\ AT 



dr 



T ' 



+ — V • F - (Per - p'Cp - pC) 



.7 



7P 

— ((7 - (Too) + KTk'^ 



AT 



(A12) 



, (A13) 



where A;^ = k"^. + k"^ is the amphtude of the total wavenumber and (Jr^ is defined by equation 
(25). Combining equations (A12) and (A13), we finally obtain equation (30) as a local 
dispersion relation for nonradial thermal perturbations. 
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Table 1. Parameters and Time Scales for Five Clusters. 



Name 


T(0) (keV)^ 


ne(0) (cm ^)'' 


r 


tcooi (Gyr)'^ 


igrow (Gyr)*^ 


/ /-/. c,d 
''grow/ ''00,0 


Abell 1795 


2 


0.049 


0.2 


0.98 


4.1 


6.3 


Abell 1835 


5 


0.17 


0.4 


0.45 


1.9 


6.3 


Abell 2199 


1.6 


0.074 


0.4 


0.58 


3.5 


9.0 


Abell 2390 


4 


0.069 


0.3 


0.98 


4.6 


7.1 


RXJ 1347.5-1145 


6 


0.11 


0.3 


0.76 


3.3 


6.2 



'^Adopted from ZN03. 

'^^cooi is the isobaric cooling time at the cluster center (see eq. [8]). 
Hgrow is the growth time of the global radial mode. 

Hoo,o is the growth time of the local isobaric thermal instabihty at the cluster center in the 
absence of conduction (see eq. [27]). 
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Fig. 1. — (a) Density, (b) temperature, (c) radial heat flux, and (d) isobaric cooling time 
(^coob solid line, see eq. [8]) of a model of A1795 with conduction. The chosen parameters 
are / = 0.2, T(0) = 2 keV, ne(0) = 0.049 cm~3, Mq = 6.6 x 10^^ M©, and r, = 20rc = 460 
kpc (ZN03). Shown also in (rf) are the growth times of local {too, dashed, see eq. [27]) and 
global (dotted) thermal instability. 
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Fig. 2. — Eigenfunctions of the radial unstable mode for five clusters, plotted as functions 
of radius. Note that the amplitudes of the solutions are largest near the center and decrease 
very rapidly as r increases. 
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Fig. 3. — Effect of thermal conductivity on the growth rates of global unstable modes in 
A1795. The abscissa is the reduction factor / of thermal conductivity applied to pertur- 
bations, while the ordinate is the eigenfrequency normalized by croo,o = (0.64 Gyr)~^, the 
growth rate of the local isobaric instability at the cluster center. Each curve is labeled by n, 
the number of nodes in the corresponding eigenfunction $,r,n/f- See text for details. 
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Fig. 4. — (a) Examples of eigenfunctions for decaying modes {a < 0) in A1795 with / = 0.2 
{b) The effect of thermal conductivity on the decay rates of modes. Higher-order modes have 
larger decay rates. 




Fig. 5. — Time evolution of (a) maximum density and (b) mass inflow rate measured at 
r = 10 kpc in models A-C. The two dotted lines in (a) represent the maximum growth rates, 
a{f = 0.2) = (4.1 Gyr)^^ and o"oo,o = (0.64 Gyr)^-*^, of thermal instability with and without 
conduction as estimated from the linear analysis. As expected, these lines agree well with the 
time evolutions of Models B and C. The insert in (a) compares the evolution of the maximum 
density in model A (sohd line) with the analytic prediction, ne(t)/ne(0) = {l — 3t/2tcooi) '^^^, 
for isobaric cooling (dotted line). See text for details. 
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Fig. 6. — Evolution of {top row) electron number density, {middle row) radial velocity in 
units of the central sound speed cq = 390 km s^^, and {bottom row) temperature, for models 
A-C. Left: Strong radiative cooling in model A quickly leads to an inward cooling flow, which 
causes the central density to increase by more than three orders of magnitude in 0.9 Gyr. 
Middle: Model B, which is initially in equilibrium but has no conduction, has perturbations 
that grow at a rate (Too, driving the system into a highly nonlinear state in less than 5 
Gyr. Right: Thermal conduction in model C delays the development of thermal instability 
significantly. Even after 20 Gyr, the density and temperature profiles are not very different 
from the initial ones. 



